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Abstract 

This paper is to introduce a type of full multigrid method for the nonlinear 
eigenvalue problem. The main idea is to transform the solution of nonlin¬ 
ear eigenvalue problem into a series of solutions of the corresponding linear 
boundary value problems on the sequence of finite element spaces and nonlin¬ 
ear eigenvalue problems on the coarsest finite element space. The linearized 
boundary value problems are solved by some multigrid iterations. Besides the 
multigrid iteration, all other efficient iteration methods for solving boundary 
value problems can serve as the linear problem solver. We will prove that the 
computational work of this new scheme is truly optimal, the same as solving 
the linear corresponding boundary value problem. In this case, this type of 
iteration scheme certainly improves the overfull efficiency of solving nonlinear 
eigenvalue problems. Some numerical experiments are presented to validate 
the efficiency of the new method. 
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1 Introduction 


In recent years, much effort has been devoted to the study of problems in solving 
large scale eigenvalue problems. Among these eigenvalue problems, there exist many 
nonlinear eigenvalue problems [3, 4, 10, 11, 12, 13, 16, 18, 20, 23], for instance the 
calculation of the Gross-Pitaevskii equation describing the ground states of Bose- 
Einstein condensates [3, 4] or the Hartree-Fock and Kohn-Sham equations used to 
calculate ground state electronic structures of molecular systems [11, 12, 16, 18, 20, 
23] from physics, chemistry and material science. However, these high-dimensional 
eigenvalue problems are always very difficult to solve. 

The multigrid and multilevel methods [2, 5, 6, 7, 8, 15, 19, 21, 22, 24, 29] provide 
optimal order algorithms for solving boundary value problems. The error bounds 
of the approximate solutions obtained from these efficient numerical algorithms are 
comparable to the theoretical bounds determined by the finite element discretization. 
But there is no many efficient numerical methods for solving nonlinear eigenvalue 
problems with optimal complexity. Recently, a type of multigrid method for eigen¬ 
value problems has been proposed in [17, 25, 26, 27, 28]. The aim of this paper is 
to present a full multigrid method (sometimes also referred as nested finite element 
method) for solving nonlinear eigenvalue problems based on the combination of the 
multilevel correction method [25, 27] and the multigrid iteration for boundary value 
problems. Comparing with the method in [17, 25, 27, 28], the difference is that it 
is not necessary to solve the linear boundary value problem exactly in each correc¬ 
tion step. We only get an approximate solution with some multigrid iteration steps. 
In this new version of multigrid method, solving nonlinear eigenvalue problem will 
not be much more difficult than the multigrid scheme for the corresponding linear 
boundary value problems. 

An outline of the paper goes as follows. In Section 2, we introduce the finite 
element method for eigenvalue problem and state some basic assumptions about 
the error estimates. A type of full multigrid algorithm for solving the nonlinear 
eigenvalue problem and the corresponding computational work estimate are given 
in Section 3. Two numerical examples are presented in section 4 to validate our 
theoretical analysis. Some concluding remarks are given in the last section. 


2 Finite element method for nonlinear eigenvalue 
problem 

This section is devoted to introducing some notation and the finite element method 
for nonlinear eigenvalue problem. In this paper, the standard notation for Sobolev 
spaces IT*’^(r2) and their associated norms and semi-norms (cf. [1]) will be used. 
For p = 2, we denote and i7o(G) = {n G : v\qq = 0}, 

where n|n = 0 is in the sense of trace, || • = || • ||s, 2 ,o. Let V = and 
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II ■ IIs denote || • ||s^n for simplicity. To facilitate the following instructions, the letter 
C (with or without subscripts) denotes a generic positive constant which may be 
different at its different occurrences through the paper. 

This paper is concerned with the following nonlinear elliptic eigenvalue problem: 
Find {X,u) G TZ X Hq{Q) such that 

{ —'V-{AVu) + f{x,u) = Xu, in fl, 

u = 0, on (2.1) 

J^u^dn = 1 , 

where ^ is a symmetric and positive dehnite matrix with suitable regularity, /(x, u) 
is a nonlinear function corresponding to the variable u, Q G TZ^ {d = 2, 3) is a 
bounded domain with Lipschitz boundary dQ. 

In order to use the hnite element method for the eigenvalue problem (2.1), we 
dehne the corresponding variational form as follows: Find {X,u) E TZ x V such that 
b{u, u) = 1 and 


a{u,v) = Xb{u,v), Wv E V, (2.2) 

where 

a{u,v) = j (^AVu ■ Vv + f{x,u)v^dil, b{u,v) = J uvdVt. 

For simplicity of describing and understanding, we only consider the numerical 
method for the simple eigenvalue case. 

Now, let us dehne the hnite element approximations for the problem (2.2). First 
we generate a shape-regular decomposition of the computing domain C 7Z^ {d = 
2,3) into triangles or rectangles for d = 2 (tetrahedrons or hexahedrons for d = 3) 
(cf. [9, 14]). The diameter of a cell K E Th is denoted by and the mesh size h 
describes the maximum diameter of all cells K E T- Based on the mesh Th, we can 
construct a hnite element space denoted by 14 C V. For simplicity, we set 14 as the 
linear hnite element space which is dehned as follows 

Vh = {vh E c{n) I Vh\K e Pi, VK e %}, (2.3) 

where Vi denotes the linear function space. 

The standard hnite element scheme for eigenvalue problem (2.2) is: Find {Xh, Uh) G 
TZxVh such that b{uh,Uh) = 1 and 


a{uh,Vh) = Xhb{uh,Vh), ^VhEVh. 


(2.4) 


Dehne a bilinear form a(-, •) as follows 


d{w, v) 


AVw ■ VvdVt, 


WweV, yvEV 
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and the correspoding norm || ■ ||a is defined by 


ll^lla = ^/a{v,v). (2.5) 

Denote 

dh{u) = inf \\u-Vh\\a- (2.6) 

Vh&Vh 

For designing and analyzing the full multigrid method, we state the following 
assumption for the nonlinear function f{x, ■) : TZ x V ^ V. 

Assumption A; The nonlinear function f{x, •) has the following estimate 
\{f{x,w) - f{x,v),ij)\ < Cf\\w - v\\o\\v\\a, e V, \/v e V, yi/j e V. (2.7) 

For generality, we only state the following assumptions about the error estimate for 
the eigenpair approximation {Xh,Uh) defined by (2.4) (see, e.g., [10, 12] for practical 
examples). 

Assumption Bl: The eigenpair approximation (A/i, Uh) of (2.4) has the following 
error estimates 


\\u-Uh\\a < {I + Cur]a(yh))Sh{u), (2.8) 

|A - Ah| + ||m - h/illo < Cur]aiVh)\\u - Uh\\a, (2.9) 

where r^aiVh) depends on the finite dimensional space 14 and has the following 
property 

= Va{Vh) < Va{Vh) iiVhCVhCV. (2.10) 

h^O 

Here and hereafter Cu is some constant depending on regularity of mesh and the 
exact eigenfunction but independent of the mesh size h. 

Assumption B2: Assume is a subspace of 14. Let us define the eigenpair 
approximation by solving the nonlinear eigenvalue problem as follows: 

Find (A^,n^) eTZxV^ such that b{u^,u^) = 1 and 

a(M^^;'^) = A'*fe(M^n'*), e (2.11) 

Then the following error estimates hold 

\\Uh-U^\\a < {^ + Cur]a{V^))bh{Uh)-: (2.12) 

|A,-A"| + ||h,-n"||o < aha(l^")||h/.-w"IU, (2.13) 

where 

5h{uh) ■= inf ||%-n''||i. (2.14) 
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3 Full multigrid algorithm for nonlinear eigen¬ 
value problem 


In this section, a type of full multigrid method is presented. In order to describe 
the full multigrid method, we hrst introduce the sequence of hnite element spaces. 
We generate a coarse mesh Th with the mesh size H and the coarse linear hnite 
element space Vh is dehned on the mesh Th- Then a sequence of triangulations 7^*, 
of n C TZ^ is determined as follows. Suppose 77i (produced from Th by regular 
rehnements) is given and let 77*. be obtained from Th^_^ via one regular rehnement 
step (produce subelements) such that 




^ hk—i, 


k = 2, - ■ ■ , n. 


(3.1) 


where the positive number /3 denotes the rehnement index and larger than 1 (always 
equals 2). Based on this sequence of meshes, the corresponding nested linear hnite 
element spaces can be built such that 


c---cW„. (3.2) 

The sequence of hnite element spaces 14^ C 142 C • • • C 14^ and the hnite element 
space Vh have the following relations of approximation accuracy (cf. [9, 14]) 

VaiVn) >C6h,{u), 6h^{u) = k = 2,---,n. (3.3) 


3.1 One correction step 


In order to design the full multigrid method, we hrst introduce an one correction 
step in this subsection. Assume we have obtained an eigenpair approximation 
^ TZ X 14^, where {£) denote the 4th iteration step in the A:-th level 
hnite element space 14^,. In this subsection, a type of correction step to improve the 
accuracy of the current eigenpair approximation ^e given as follows. 

Algorithm 3.1. One Correction Step 

1. Define the following auxiliary boundary value problem: Finduf^^'^ G 14^, such 
that 


a[u 






(3.4) 


eigenfunction approximation ^ 


/£'\ 

Perform m multigrid iteration steps with the initial value to obtain a new 


u 


(^+ 1 ) 


= MCiVk^Xhl^hl - 


h}^ / 




(3.5) 
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(i) (i) 

where Vh^. denotes the working space for the multigrid iteration, — 

f{x,ufl_) is the right hand side term of the linear eguation, denotes the 
initial guess and m is the number of multigrid iteration times. 

2. Define a new finite element space Vn^hk = and solve the 

following eigenvalue problem: Find E TZ x VH,hk such that 


b{u 


(£+1) „,(£+!) 




u 




) = 1 and 


‘'hk ’ 


= A 


^k 


'^VH,hk e (3.6) 


In order to simplify the notation and summarize the above two steps, we define 


(A, 


(^+1) 

^k 




= EigenMGiVn, A^, Vh^,m) 


Theorem 3.1. Assume the multigrid iteration m) of 

( 3 . 4 ) has the following error reduction rate 


\u 


.(£+ 1 ) 


— u 


(^+ 1)1 

hk I 


a < “W^hk lU’ 


.W| 


and the given eigenpair approximation has following estimates 


I A,.. - A'" 




+ - “£llo < C«'(a(VH)|ItH, - 


(3.7) 


(3,8) 


Under Assumptions A and B2, the resultant eigenpair approximation G 

IZ X V/ij, produced by performing Algorithm 3.1 has the following error estimates 

-“I'A’llo ^ - Al'llo. (3,9) 

|A». - Ai'^’l + < CMVH)\\aH, - (3.10) 

where 

'^ = 0 F (^CuO + (1 + (^Cu + C/) (1 + CuPa{VH))'j PaiVn) (3-11) 

and Cu depends on the desired eigenpair. 


Proof. From (2.4), (2.7) and (3.4), we have 

a{Uhk - = (^{XhkUhk - Uhf) - fix, M^)) , Vh,^ , 

< \K\\\uhk - UhlWoWvhkWo + \K - ^hlWWhlW 
+Cf\\hhk - Uhl\\o\\vhk\\a, Vn/,, e Vh^. 

It leads to the following estimates by using the property of a(-, •) and (3.8) 

\\Uhk-U^h^^^\\a < iCu + Cf)paiVH)\\Uh^-U^hl\\a, (3-12) 
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where Cu depends on the desired eigenpair. 

Combining (3.7) and (3.12) leads to the following error estimate for 


I ^h. l a 








(t+1) 


u 


W|l 'i 
hk 11 “/' 

W| 


— ^{W^hk '^hkWa \\Uhf, 

< 9{l + {C^ + Cf)va{VH))\\u,,-u\:i 
Then from (3.12) and (3.13), we have the following inequalities 


(3.13) 


I- ~h+i)i 


< IL-;. _JT(^+b| 


S Uh,-U 




+ iitth+i) _w(^+i)| 
^ W^hk “/ifc I 


< (0 + (1 + e){Cu + CMVh))\\ uh, - uflWa. (3.14) 


The eigenvalue problem (3.6) can be regarded as a finite dimensional subspace 
approximation of the eigenvalue problem (2.4). Using (2.12) and (2.13) in Assump¬ 
tion B2, the following estimates hold 


\Uh^ - Uh 


h+i)i 


< {l + Cur]a{VH,hk)) inf \\Uh;, - VH,hk\\a 

< (1 + Cur]a{VH)) 


(^+ 1 ) I 


< l\\Uh,-ufl\\a, 


(3.15) 


and 


|A/.,-A 


h+i) I 
I 


+ - “!.7'’iio £ c„i7„(Vh,).jiimih - 47'’' 


(3.16) 


Then we obtained the desired results (3.9) and (3.10) and the proof is complete. □ 


3.2 Full multigrid method for eigenvalue problem 

In this subsection, based on the one correction step defined in Algorithm 3.1, a type 
of full multigrid scheme will be introduced. The optimal error estimate with the 
optimal computational work will be deduced for this type of full multigrid method. 

Since the multigrid method for the boundary value problem has the uniform error 
reduction rate (cf. [9, 15]), we can choose suitable m such that 0 < 1 in (3.7). From 
the definition (3.11) for 7 , it is obvious that 7 < 1 when the mesh size H of Th is 
small enough. Based on these property, we can design a full multigrid method for 
nonlinear eigenvalue problem as follows. 

Algorithm 3.2. Full Multigrid Scheme 
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1. Solve the following nonlinear eigenvalue problem in Vh^: Find {Xhi,Uhi) G 

X V/jj such that b{uhi,Uhi) = 1 and 

'^hf) Xfi^bl^Ufi^ , Vh^f V'y/jj G Vhi' 

Solve this nonlinear eigenvalue problem to get the desired eigenpair approxi¬ 
mation {Xh^.Uhf) eTZx V/ij. 

2. For k = 2, - ■ ■ ,n, do the following iterations 

• Set Xi°J = Xh,_, and = Uh,_,. 

• Perform the following multigrid iterations 

= EigenMCiVn, A^, K,, m), for £ = 0, • • • ,p - 1. 

• Set Xh, = AifJ and Uh, = 4^- 
End Do 

Finally, we obtain an eigenpair approximation {Xh„,Uh,,) E 71 x Vh„ in the finest 
space. 

Theorem 3.2. Assume the conditions of Theorem 3.1 and Assumption B1 hold. 
After implementing Algorithm 3.2, the resultant eigenpair approximation {Xh„,Uh„) 
has the following error estimate 

\\uh„-UhAa < (3.17) 

IV - AftJ + \\uh^ - Uhjo < ^^^p VaiVH)5hM, (3.18) 

under the condition (3^^ < 1 . 


Proof. Define := Uhf, — Uh^.. Then from step 1 in Algorithm 3.2, it is obvious 
Cl = 0. For k = 2, - ■ ■ ,n, from Assumption B1 and Theorem 3.1, we have 

||cA;||a ^ ^/ifc_i||a 

< 7^{\\Uh, - Uh,_,\\a + \\Uh,_, - Uh,_,\\a) 

< 'yP{C6h^{u)+ \\ek-i\\a)- (3.19) 

By iterating inequality (3.19) and the condition < 1, the following inequalities 
hold 


|e„|U < C7^4„(n) + C72^4„_,(n) + --- +V-'V,(n) 

n / ^ 


, n—k 


k =2 


, fc =2 


< 


7^4 J«) 




(3.20) 


< C- 


T 


1 - /37P 


- 5 , 


[u] 


For such choice of p, we arrive the desired result (3.17) and (3.18) can be obtained 
by (2.13), (3.10) and (3.17). □ 


Remark 3.1. The good convergence rate of the multigrid method for boundary value 
problems leads to that we do not need to choose large m and p [9, 15, 22, 29]. 


3.3 Estimate of the computational work 

In this subsection, we turn our attention to the estimate of computational work for 
the full multigrid method defined in Algorithm 3.2. It will be shown that the full 
multigrid method makes solving the nonlinear eigenvalue problem need almost the 
same work as solving the corresponding linear boundary value problems. 

First, we define the dimension of each level finite element space as '■= diml4j.. 
Then we have 

/ 1 \ d{n-k) 

iVn, k = l,2,---,n. (3.21) 

The computational work for the second step in Algorithm 3.1 is different from 
the linear eigenvalue problems [17, 25, 26, 27]. In this step, we need to solve a non¬ 
linear eigenvalue problem (3.6). Always, some type of nonlinear iteration method 
(self-consistent iteration or Newton type iteration) is adopted to solve this nonlinear 
eigenvalue problem. In each nonlinear iteration step, it is required to assemble the 
matrix on the finite element space {k = 2, - ■■ ,n) which needs the computa¬ 

tional work 0{Nk). Fortunately, the matrix assembling can be carried out by the 
parallel way easily in the finite element space since it has no data transfer. 

Theorem 3.3. Assume we use d computing-nodes in Algorithm 3.2, the nonlin¬ 
ear eigenvalue solving in the coarse spaces VH,hk (k = 1, • • • > W ^1 need work 
0 {Mh) and 0{Mhf), respectively, and the work of the multigrid solver MGfVh,., 

m) in each level space 14^, is 0{Nk) for k = 2,3, ■ ■ ■ ,n. Let vj denote the non¬ 
linear iteration times when we solve the nonlinear eigenvalue problem (3.6). Then 
in each computational node, the work involved in Algorithm 3.2 has the following 
estimate 

Total work = 0^(l + ^)Nn Mh log A'n + ■ (3.22) 

Proof. We use M4 to denote the work involved in each correction step on the k-th 
finite element space 14^.. From the dehnition of Algorithm 3.1, we have the following 
estimate 

Wk = 0 fNk + Mj, + (3.23) 
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Based on the property (3.21), iterating (3.23) leads to 


Total work = ^ ^ hhfc — Cl ( Tf/ij + ^ ^ f 


k=l 


k=2 




c? I (l + + (ri - 1)Mh + M,, 


,/c=2 


2 \ d(n—A:) 


O ((1 + j)N^ + MnlogN^ + M;,,) . 


This is the desired result and we complete the proof. 


(3.24) 

□ 


Remark 3.2. Since we have a good enough initial solution Uhf.j^^ in the second step 
of Algorithm 3.1, then solving the nonlinear eigenvalue problem (3.6) always does 
not need many nonlinear iteration times (always zu < 3). In this case, the complexity 
in each computational node will be 0{Nn) provided Mh <C Nn and < N^. 


4 Numerical results 

In this section, two numerical examples are presented to illustrate the efficiency of 
the full multigrid scheme proposed in this paper. 

Example 4.1. In this example, we consider the ground state solution of Gross- 
Pitaevskii eguation (GPE) for Bose-Einstein condensation (BEG). 

{ —Au + Wu + (\u\‘^u = Xu, in n, 

u = 0, on dll, (4.1) 

J^u^dfl = 1, 

where hi denotes the three dimensional domain [0,1]^ , C = 1 o.nd W = x1 + X2 + x^. 

From the results [10, 28], the Assumptions A, B1 and B2 hold for the GPE 
(4.1). So the proposed full multigrid method can be applied to the GPE (4.1). The 
sequence of hnite elements spaces are constructed by linear element on a series of 
meshes produced by regular refinement with /3 = 2. In each level of the full multigrid 
scheme defined in Algorithm 3.2, the parameters are set to be m = 1, p = 1. And 
we take 3 conjugate gradient smooth steps for the presmoothing and postsmoothing 
iteration step in the multigrid iteration in the step 1 of Algorithm 3.1. Since the 
exact solution is not known, an adequate accurate approximation is chosen as the 
exact solution for our numerical test. Figure 1 shows the corresponding initial mesh. 
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Figure 1: The initial mesh for Examples 4.1 

Figure 2 gives the corresponding numerical results of Algorithm 3.2. From Figure 
2 , we can hud that the full multigrid scheme can obtain the optimal error estimates 
for both eigenvalue and eigenvector. 

In order to show the efficiency of Algorithm 3.2, we provide the CPU time for 
Algorithm 3.2. Here, we choose the Package ARPACK as the eigenvalne solving tool 
and the full multigrid scheme is running on the machine PowerEdge R720 with the 
linux system. The corresponding results are presented in Table 1 which shows the 
efficiency and linear complexity of Algorithm 3.2. 


Table 1: The CPU time for Example 4.1 by Algorithm 3.2 


Number of levels 

Number of elements 

Time for Algorithm 3.2 

1 

3072 

0.45 

2 

24576 

1.55 

3 

196608 

8.08 

4 

1572846 

63.01 

5 

12582912 

519.86 


Example 4.2. In the second example, we consider the GPE with the coefficient 
( = 100 and W = xl + xl + on the domain = [0,1]^. 

The initial mesh used in this example is the one shown in Figure 1. Numerical 
results are present in Table 2 and Figure 3. It is obvious. Table 2 and Figure 3 also 
show that the efficiency and linear complexity of Algorithm 3.2. 


5 Concluding remarks 

In this paper, a type of full multigrid method is introduced for nonlinear eigenvalue 
problems. The proposed methods is based on the combination of the multilevel 
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Errors by full multigrid method 


CPU time of full multigrid method 




Figure 2: Left; The errors of the full multigrid method for the ground state solution 
of GPE, where and Uh denote the numerical eigenvalue and eigenfunction by 
Algorithm 3.2. Right: CPU Time of Algorithm 3.2 for Example 4.1 


Table 2: The CPU time for Example 4.2 by Algorithm 3.2 


Number of levels 

Number of elements 

time for Algorithm 3.2 

1 

24576 

4.32 

2 

196608 

11.43 

3 

1572846 

70.88 

4 

12582912 

577.52 


correction technique for nonlinear eigenvalue problems and the multigrid iteration for 
linear boundary value problems. The multilevel correction technique can transform 
the nonlinear eigenvalue solving into a series of solutions of linear boundary value 
problems on a sequences of hnite element spaces. The multigrid iteration is one of 
the efficient iteration which has uniform error reduction rate. 

The multigrid iteration can also be replaced by other types of efficient iteration 
schemes such as algebraic multigrid method, the type of preconditioned schemes 
based on the subspace decomposition and subspace corrections (see, e.g., [9, 29]) 
and the domain decomposition method (see, e.g., [24, 30]). 
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